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The inhomogenous real-space electronic structure of gapless and gapped disordered bilayer 
graphene is calculated in the presence of quenched charge impurities. For gapped bilayer graphene 
we find that for current experimental conditions the amplitude of the fluctuations of the screened 
disorder potential is of the order of (or often larger than) the intrinsic gap, A, induced by the 
application of a perpendicular electric field. We calculate the crossover chemical potential, Acr, 
separating the insulating regime from a percolative regime in which less than half of the area of 
the bilayer graphene sample is insulating. We find that most of the current experiments are in the 
percolative regime with Acr ^ A. The huge suppression of Acr compared with A provides a possi- 
ble explanation for the large difference between the theoretical band gap A and the experimentally 
extracted transport gap. 

PACS numbers: 



One of the unique properties of single layer graphene, 
SLG, [1] is its high room-temperature electronic mobility 
[2]. This fact makes graphene of great interest for pos- 
sible technological applications. However, the lack of a 
band gap implies that in SLG the current can never be 
turned off completely (i.e. SLG has a very low on/off 
ratio in the engineering jargon) and therefore limits the 
possible use of SLG in transistor or switching applica- 
tions. The mobility of bilayer graphene, BLG, is nor- 
mally lower than the mobility of SLG [IfE], but it can 
be very high when Boron-Nitride is used as a substrate 
[6l [7 . Most importantly, by applying a perpendicular 
electric field P4T4] , a gap of up to 250 meV can be opened 
in the band structure of BLG which should strongly en- 
hance the on/off switching ratio. The strictly 2D na- 
ture of the carriers, the high room-temperature mobil- 
ity, and the ability to open and tune A, make BLG an 
extremely interesting material both from a fundamental 
and a technological point of view. In recent BLG ex- 
periments pTHT3] the activated transport gap has been 
found to be orders of magnitude smaller than A. This 
finding has been explained assuming that transport is in 
the variable range hopping regime, VRH, pTtiT3] . or that 
edge modes might contribute significantly to transport 
[15]. Resonant scattering centers have been proposed as 
the dominant source of disorder [16 in both SLG and 
BLG. In gapped BLG the resonant scatterers would in- 
duce localized states that would then mediate transport 
via VRH. However scanning tunneling microscopy exper- 
iments have so far not shown direct evidence of resonant 
states suggesting that their density might be quite low, 
in addition no sign of localization is ever observed in un- 
gapped BLG even in the presence of strong disorder. On 
the other hand, the edge modes can significantly con- 
tribute only if scattering between counter-propagating 
edge modes is suppressed. Because the BLG band struc- 



ture is characterized by two equivalent valleys with op- 
posite chirality, even small quantities of short-range de- 
fects can mix the fermionic states of the two valleys and 
greatly suppress the contribution of the edge modes to 
transport [14 . These facts motivated us to look for a 
possible alternative explanation to the smallness of the 
experimental BLG transport gap compared to A based 
on the disorder-induced massive breakdown of momen- 
tum conservation. 

Our model is based on the assumption that charge im- 
purities are the dominant source of disorder in exfoliated 
BLG samples. There is ample evidence [2 that this as- 
sumption is at least consistent with most of the transport 
experiments on gapless SLG and BLG although other 
scattering sources might play an important role [16l [17] . 
In particular we assume the charge impurities to be lo- 
cated at a typical distance d ^ 1 nm from the graphenic 
layer and to be uncorrelated. In reality some degree of 
correlation is expected but it does not affect qualitatively 
our results. In the presence of charge impurities the car- 
rier density becomes strongly inhomogenous. The impor- 
tance of density inhomogeneities for the understanding 
of the physics, especially transport, of 2D electronic sys- 
tems has been appreciated in the context of the Quan- 
tum Hall, QH, effect [18 and of standard 2D electron 
gases, 2DEGs, in which the gap between the hole-band 
and the electron-band is much larger than the disorder 
strength The case of gapped BLG is different from 
these cases because no magnetic field is present and so 
the dispersion is not broken-up in Landau Levels and yet 
the band-gap is small compared to the strength of the 
disorder potential. We emphasize, moreover, that the 
transport phase diagram, which is the main topic of our 
work , has never before been addressed in the literature 
for any system. In BLG we therefore have the unique 
condition of a small band gap between non-degenerate 
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valence and conduction bands. The most striking and 
counterintuitive consequence of this fact is that, as we 
show below, in gapped BLG an increase of the disorder 
strength can drive the system from being an insulator to 
be a bad metal when the chemical potential, /i, is within 
the gap. This behavior is the opposite of what happens 
in a standard 2DEG in which an increase of the disor- 
der strength drives the system to an insulating state. By 
showing that the disorder can effectively drive BLG to a 
metal state even when A is finite and ji is well within the 
gap, our work provides a compelling possible explanation 
for the large discrepancy between A of gapped BLG and 
the gap extracted from transport measurements. We de- 
fine an effective real -space gap, A^r^ that determines the 
transport properties and show that for disorder strengths 
typical in current experiments it can be zero even for A 
as large as 150 meV. Although our specific calculations 
are carried out for gapped BLG systems, the general idea 
developed here, namely a spatially fiuctuating local band 
gap in the presence of charged impurity disorder induced 
inhomogeneity, should apply to other systems, and we 
believe that the same idea could explain the experimen- 
tal finding [20 of rather small transport band gaps in 
graphene nanoribbon experiments where percolation ef- 
fects are known to be important. 

Our main goal is to find a qualitative explanation for 
the smallness of the transport gap compared to A for the 
situation when charge impurities are the dominant source 
of disorder. To calculate the electronic structure in the 
presence of the disorder potential due to charge impu- 
rities we use the Thomas-Fermi theory, TFT. For SLG 
the TFT results [21] compare well with density functional 
theory, DFT, results [22] and experiments [23H25] as long 
as the impurity density nimp is not too low lO^^cm"^) 
[26] . These results suggest that TFT might give rea- 
sonable results also for disordered BLG. TFT is valid 
when the density profile, n(r), satisfies the inequality 
\Vn/n\ < kp^ with kp = ^/7Tn the Fermi wavevector. 
As shown below the density varies on length scales of 
the order of 10-20 nm whereas in the metallic regions 
n ^ 2 X lO^^cm"^, so that the inequality |Vn/n| < /c^? is 
only marginally satisfied. A complete quantitative vali- 
dation of the TFT results can only be achieved by com- 
parison to DFT results that however are not yet available 
for disordered BLG. On the other hand due to the com- 
putational cost DFT cannot be used to calculate disorder 
averaged quantities that are needed to extract the trans- 
port properties. TFT is therefore the only approach that 
can be used to address the issue of transport in BLG 
in the presence of long range disorder and, given that 
our goal is the qualitative understanding of the large dif- 
ference between transport gap and A, is also adequate. 
Moreover, the strong dependence of the transport gap 
on the details of the experiments (like the temperature 
range) makes a quantitative comparison of theory and ex- 
periments almost impossible. For this reason, and given 



the limited quantitative accuracy of TFT, for the band 
structure of BLG we use the simple model of 2 parabolic 
bands with effective mass m* = O.OSSmg. The simplic- 
ity of this model allows us to identify the few parame- 
ters that affect the qualitative features of the results and 
makes our findings relevant also to standard parabolic 
2DEGs. 

The TFT energy functional is given by: 
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where e the static background dielectric constant, and 
eVo/e the bare disorder potential. The first term is 
the kinetic energy, the second the Hartree part of the 
Coulomb interaction, and the third is the contribution 
due to the disorder potential. Assuming e = 4 for the 
substrate, as appropriate for graphene on Si02, in the 
remainder we set e = 2.5, average of the dielectric con- 
stant of vacuum and substrate. By differentiating E[n] 
with respect to n we find: 
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where we have introduced the screening length r^c = 
[(2e2m*)/(en2)]-i 2 nm. 

At the energy minimum, SE/5n = 0. For the gapless 
case this equation can be solved analytically in momen- 
tum space to find 
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Let Vz),sc(r) = r-'Voir) + (l/2rsc) / dr^n{r^)/\r - be 
the screened disorder potential. Using Eq. (|3| we find: 
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We have Vd{c[) = A{q)e ^^/q, where A{q) are random 
numbers with Gaussian distribution such that (A) = 
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3, where the angle brackets denote av- 



erage over disorder realizations. Using Eq. (|3|, Q and 
the statistical properties of A(q), for the variance of the 
density, Var(n), we find 
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dimensionless function (here r(a,x) is the incomplete 
gamma function). For small d/r^c^ f = —1— 7 — 
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FIG. 1: (Color online). Color plot of n(r) at the CNP 
for a single disorder realization in gapped BLG, for riimp = 
8 X 10^^ cm"^, d = 1 nm, and A = 125 meV (a) and 
A = 250 meV (b). The white areas in (a) and (b) repre- 
sent insulating regions, n^^™^-* (c) and V^^^"^ (d) as a func- 
tion of doping for riimp = 3 x 10^^ cm~^ and different val- 
ues of A, from top to bottom in (c) (bottom to top in (d)) 
A - 12,23,35,46,55 meV. The dashed line in (c) and (d) 
shows the value of n*^^™^-* and V^^^^ respectively for the gap- 
less case. The insets show the results for A = 500 meV [29]. 



log(2d/rsc) + 0{d/rsc) (where 7 = 0.577216 is the Eu- 
ler constant), whereas for d ^ r^^ f = l/(2<i/rsc)^ + 
0{{d/r^c)~^)' For the root mean square of the den- 
sity, n^™^\ and screened disorder potential, V^™^ = 



{h^ /em*)VD,sci we then have: 



(rms) 



nl/2 



(rms) 



re 



TT 



2m* 



(rms) 



(6) 



In the presence of a band gap the equation SE/Sn = 
becomes nonlinear and an analytic expression for n^^"^®^ 
is not readily obtainable. We have solved the problem 
numerically for 100 x 100 nm samples, with 1 nm spatial 
discretization, considering several (1000 or more) disor- 
der realizations to then calculate the disorder averaged 
quantities. In Fig. [l] (a) and (b) we show our calcu- 
lated carrier density landscapes for single disorder real- 
izations for two values of A, both with the same nimp 
and /i fixed at the charge neutrality point, CNP. In the 
absence of disorder both situations in Fig. [l] (a) and 
(b) will manifest zero carrier density throughout with 
a pure intrinsic band gap at all spatial points (i.e. both 
sets of plots will be completely 'white' in color since we 
are explicitly at T=0 with no thermal inter-band ex- 
citations). Fig. [1] (c) and (d) show the dependence of 
^(rms) V^^^ , respectively, on the doping n for val- 



ues of the gap between 12 and 55 meV. n*^^"^®^ is sup- 
pressed close to the CNP, due to the large area covered 
by insulating regions (n = 0) whereas V^^^^ is higher 
close to the CNP due to the lack of screening. The re- 
sults of Fig. [1] (c) and (d) are well fitted by the scalings 
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2.31e~^/^^"^^ and n and A expressed in units of 
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lO^^cm"^ and meV respectively. 

The analysis of the disorder averaged results obtained 
from the TFT can be used to qualitatively understand 
the electronic transport in the highly inhomogeneous 
density landscape of gapped BLG. Let Ai be the dis- 
order averaged fraction of the area occupied by the insu- 
lating regions (white in Figs, [l] (a) and (b)). The color 
plots in Fig. [2] (a)- (c) show the dependence of Ai on nimp 
and the doping (n) for three different values of A. At 
high dopings and low disorder Ai is fairly small. As (n) 
decreases some hole puddles and insulating regions start 
appearing. The black line identifies the contour for which 
^(rms) _ ^^^^ Below the black line n*^^"^^^ > (n) and we 
are in the strongly inhomogenous regime. The white line 
identifies the contour for which the area fraction occupied 
by electron puddles is equal to 50%. Below this line more 
than half of the sample area is occupied by hole puddles 
and insulating regions. Finally, the red line shows the 
contour for which is Ai = 0.5. Below this line more than 
half of the sample is covered by insulating regions. The 
counterintuitive result is that close to the CNP Ai be- 
comes smaller as nimp increases with (n) fixed. This is a 
consequence of the fact that as nimp increases, the disor- 
der becomes strong enough to bring the Fermi level inside 
the conduction or the valence band. Using the contour 
lines overlaid on the color-plot for Ai we can qualitatively 
identify different semiclassical transport regimes. Regime 
I: Weak disorder where n^^"^®^ <C (n). The system is a 
good metal with an almost uniform density landscape. 
Regime II: Strong disorder n'^^^^^ > (n). In this regime 
the density landscape breaks up in puddles and insulating 
regions. Because the area fraction covered by electrons is 
larger than 1/2 the system behaves like a metal via per- 
colation through the electron puddles. Regime III: Weak- 
moderate disorder, < A/2, 1/^™^ < A/2. Most of 
the sample area is insulating and the system is an insula- 
tor. Regime IV: In this regime the disorder is so strong, 
y^™^ ^ ^ ImI' that, despite the finite gap, Ai is 
less than 50% but neither the electron puddles, nor the 
hole puddles alone cover more than 50% of the total area. 
The system should behave as a bad metal with the con- 
ductance determined by tunneling events across few nar- 
row insulating regions separating the electron-hole pud- 
dles similar to transport at the CNP in gapless BLG [27j 
and SLG [28]. It is obvious, looking at the color plots in 
Fig. 2(a)-(c) that at lower (higher) values of A (disorder) 
the apparent transport gap will be strongly suppressed 
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FIG. 2: (Color online). Color plots showing Ai as a func- 
tion of riimp and (n) for A = 125 meV, (a), A = 250 meV, 
(b), and A = 500 meV, (c). (d) Acr as a function of 
rzimp for different values of A. (e) Acr as function of A 
for different values of riimp- (f) Acr/(A/2) as a function of 
Vd/A. The filled (unfilled) symbols, connected by the solid 
(dashed) line, show the dependence of Acr/(A/2) with re- 
spect to V^'^;^/A ((ey^™'Ve)/A). The different symbols, 
circles, squares, diamonds, triangles, show the results for 
A = 500, 250, 125, 62.5 meV respectively. 29 . 



and in fact will be vanishingly small. We believe that this 
is the current physical situation in existing BLG systems. 

To identify the transport regime is useful to intro- 
duce the critical value of chemical potential, Acr, for 
which ^4^=50%. For /i > A^ (/i < Acr) the system 
is expected to behave as bad metal (an insulator). For 
A Acr ^ A/2 is almost independent of the dis- 

order strength. For small gaps, A < Acr depends 

strongly on the disorder strength, and Acr <C A. This is 
shown in Fig.[2](d) in which the calculated Acr as a func- 



tion of nimp is shown for four different values of A. For 
fixed A, Acr decreases as nimp increases. We see that for 
A < 150 meV, for impurity densities of the order of the 
ones estimated in current experiments on exfoliated BLG 
(nimp ^ 10^^ cm~^), Acr can be orders of magnitude 
smaller than A. Fig. [5] (e) shows Acr as a function of A 
for different values of nimp- We see that for A ^ 100 meV 
Acr ~ already for nimp ^ 3 lO^^cm"^. In the presence 
of spatial correlations among impurities the density in- 



homogeneities are expected to be reduced and therefore 
the value of nimp for which Acr would increase. Fi- 
nally Fig. [2] (f) shows the scaling of Acr on the strength 
of the disorder potential, Vd- The points connected by 
solid (dashed) lines show the dependence of Acr/(A/2) 
with respect to F^^^'Va {{eV^'"^'^ / e) / A) , with V^'"^'^ 
the rms of the bare disorder potential. In both cases, by 
normalizing both Acr and Vd with the band-gap, we find 
that the results obtained for different values of A col- 
lapse on a single curve that does not depend on A and 
away from Vd = scales approximately as 1 — ae^^^/^, 
with (a = 0.02,6 = 1.52) for Vd = eV^^V^, and 
(a = 0.015, b = 6.65) for Vd = V^'^'\ 

In summary, using a simple 2-bands model and TFT, 
we have characterized the density inhomogeneities in 
gapless and gapped BLG. For gapless BLG we have 
found analytic expressions for n^^"^^) and V^™^ as a 
function of the experimental parameters and shown 
that they do not depend on the doping and scale like 
V^imp/'/'(^Asc)/rsc. For gapped BLG n^^-^) (F^^^^) 
is reduced (enhanced) with respect to the gapless case, 
in particular in the vicinity of the CNP. By calculating 
the disorder averaged fraction of the sample area, A^, 
covered by insulating regions we have qualitatively iden- 
tified four different transport regimes. We have shown 
that most of the current experiments are expected to be 
in a regime, regime /F, in which the disorder is strong 
enough to reduce Ai below 50% even at zero doping. In 
this regime, gapped BLG is expected to behave like a 
bad metal in which transport is dominated by hopping 
processes between electron and hole puddles that cover 
most of the sample. The value of the chemical potential 
Acr for which Ai = 0.5 identifies the crossover between 
the insulating regime and regime IV. We have shown how 
Acr depends on the impurity density nimp, the strength of 
the screened disorder potential, and the theoretical band 
gap A. We believe the reduction of Acr as a function 
of 

^imp is the qualitative resolution of the contradiction 
between A and the transport gap. A clear prediction of 
our theory is that in cleaner BLG (e.g. on boron nitride 
substrates) there should be close agreement between the 
transport gap and A. We also predict agreement between 
A and the transport gap for very large values of A. 
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